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ABSTRACT 



o 

OO 

CD i The apparent dependence of detection frequency of extrasolar planets on the metal- 

hcity of their host stars is investigated with Monte Carlo simulations using a determin- 
istic core-accretion planet formation model. According to this model, gas giants formed 
Q^' and acquired their mass Mp through planetesimal coagulation followed by the emer- 

Q . gence of cores onto which gas is accreted. These protoplanets migrate and attain their 

asymptotic semi-major axis a through their tidal interaction with their nascent disk. 
j§ ! Based on the observed properties of protostellar disks, we generate Mp-a distribution. 

Our results reproduce the observed lack of planets with intermediate mass Mp = 10- 
; lOOMe and a < 3AU and with large mass Mp > lO^M® and a < 0.2AU. Based on the 

' simulated Mp-a distributions, we also evaluate the metallicity dependence of fraction of 

^ stars harboring planets that are detectable with current radial velocity survey. If pro- 

tostellar disks attain the same fraction of heavy elements which are contained in their 
host stars, the detection probability around metal-rich stars would be greatly enhanced 
because protoplanetary cores formed in them can grow to several Earth masses prior to 
their depletion. These large masses are required for the cores to initiate rapid gas ac- 
cretion and to transform into giant planets. The theoretically extrapolated metallicity 
dependence is consistent with the observation. This correlation does not arise naturally 
in the gravitational-instability scenario. We also suggest other metallicity dependence 
of the planet distributions that can be tested by on-going observations. 



Subject headings: planetary systems: formation - solar system: formation - stars: 
statics 
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1. Introduction 

One of the most important observed characteristics of known extrasolar planets is that their 
detection frequency r/j increases with the metaUicity of their host stars (Fischer & Valenti 
2003; Santos et al. 2004). Although could be significantly changed by the accretion of plan- 
ets/planctcsimals onto the host stars after they have evolved onto the main sequence, e.g., (Sandquist 
et al. 1998; Murray &: Chaboycr 2002), the lack of any dispersion among members of stellar clus- 
ters (Quillen 2002; Wilden et al. 2002) indicate that the impact of this process is very limited. The 
r/j-Z* correlation may also be interpreted as evidence that the formation probability of gas giant 
planets is greatly enhanced in metal-rich protostellar disks, in accordance with the core-accretion 
scenario for formation of giant planets, e.g., (Mizuno 1980; Bodenheimer h Pollack 1986; Pollack 
et al. 1996; Ikoma et al. 2000). 

In the conventional core-accretion scenario, heavy elements (metals) in a protoplanetary disk 
condense into grains which form km-sized planetesimals. ThroTigli gravitational interaction and 
coagulation, planetesimals evolve into cores through runaway (Grccnbcrg et al. 1978; Wcthcrill &: 
Stewart 1989; Aarseth et al. 1993; Kokubo & Ida 1996) and oligarchic growth (Kokubo & Ida 1998, 
2002). The cores' eccentricity is effectively damped by their tidal interaction with the ambient disk 
gas (Artymowicz 1993; Ward 1993). After the cores have swept up all residual planetesimals within 
a "feeding-zone width" ~ 10 Hill's radii (Kokubo h Ida 1998, 2002), their growth in the inner 
regions of the disks is stalled with an "isolation" mass Mc^iso- If a core attains the critical mass, 
acc ~ several Mq), before the disk gas depletion, rapid gas accretion onto it is initiated and it 
transforms into a gas giant planet. Since Mc^iso increases with the surface density (Sd) of dusts 
(composed of metals) in protoplanetary disks, formation of gas giants tends to be more prolific in 
a metal-rich environment. 

Here, using the deterministic planet formation model developed by Ida &: Lin (2004) (hereafter 
referred to as Paper I), we calculate the metallicity dependence of fraction (77J) of stars harboring 
planets that are detectable with the current radial velocity surveys. In §2, we briefly describe our 
model. In §3, we simulate the mass (Afp) - semimajor axis (a) distributions of extrasolar planets 
and their metallicity dependence through a series of Monte Carlo simulations. The theoretically 
constructed distributions are consistent with observed ones. Using the simulated distributions, we 
evaluate rjj. Our model explains the observed metallicity dependence of rjj. §4 is summary and 
discussions. 

2. Model 

We carry out Monte Carlo simulations to reproduce the observed Mp-a distribution of ex- 
trasolar planets. We basically follow the methods of Paper I except slightly different choice of 
parameters. Initial conditions of the Monte Carlo simulations are the a of planets and the dust and 
gas surface density (Eg and S^) at a in their nascent disks. For a given set of Sg(a) and T,^{a), 
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we numerically compute the entire formation and migration sequence of individual cores and gas 
giant planets. Here we briefly summarize our prescription. For details, see Paper I. 

Since their observational data over the relevant length scales (a few AU's) for planet-forming 
regions is not available, we scale their quantitative values with global factors /d and /g, to those of 
the empirical minimum mass nebula (MMN) model for our Solar system as 

= 1077ice/d(a/lAU)-3/2 [gcm-2], 

= 2.4 X 103/g(a/lAU)-3/2 [g cm'^], ^ ' 

where the compositional scaling parameter rjice = 1 inside the ice condensation radius [oice = 
2.7(M*/M0)2 AU and 4.2 outside it (Hayashi 1981). [Note that r/ice outside oice can be slightly 
smaller {r^ 3.0) (Pollack ct al. 1994).] The mass of the host star M^, is scaled with that of the Sun 
Mq. The dependence on the power-law index for and Sg are examined in a separate paper. To 
represent the decline in the observational signatures of protostellar disks on timescales r^ep ~ 1-10 
Myrs (Haisch et al. 2001), we assume exponential decay of Sg as 

/g = /g,oexp(-t/Tdep), (2) 

due to both viscous evolution and photocvaporation. We neglect gas replenishment onto the disk 
and any time variation in /j. We discuss the distributions of and /g.o below. 

With prescription (1), the cores' mass at time t after the formation of the disk is deduced to 
be (Paper I) 

(iiife)^.^/^/^? (i^)"" (w)""'° ikT i^)) 

(3) 

The dependence of Mc{t) on the mass m of the field planetesimals is due to the aerodynamical 
gas drag which regulates their velocity dispersion. Since e-folding growth time scale is longer for 
larger Mp, Mc{t) does not depend on initial mass of the core. In the nominal case, Eq. (3) with 
m = 10^2 g is used. If m is smaller or if radial migration of cores due to their tidal interaction with 
gas (Goldreich & Tremaine 1980; Ward 1986) and/or planetesimal disks (Murray et al. 1998; Ida et 
al. 2000) occurs, core accretion may be accelerated (Inaba et al. 2003; Tanaka & Ida 1999; Rice &; 
Armitage 2003). We also carry out the calculations with 3 times larger /smaller dMc/dt. Note that 
we do not consider the change in a of the cores due to the radial migration, because the migration 
can be a random walk (Nelson Sz Papaloizou 2004). 

Prior to severe depletion, the cores' growth is stalled when they attain an isolation mass 
(Kokubo & Ida 1998; Kominami & Ida 2002) which is given by 

M.„ = 0.1.,-/r(j|^)"'(^)""M. ,4) 

In our simulation, growth beyond the isolation after significant gas depletion is taken into account 
(Paper I). If radial migration of the cores occurs, their growth barrier can be bypassed (Tanaka 
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k, Ida 1999; Rice & Armitage 2003) up to Mc^noiso which is given with Mc^noiso — TrEda^ (feeding 
zone width Aa ~ a) by (Paper I) 



In the nominal case, the asymptotic masses of the cores are given by min(Mc,iso, Mc^noiso)- We also 
examine tha case where Mc,noiso is adopted in entire /d and a. 

For sufficiently rapid gas accretion to transform a core into a gas giant planet, it must attain a 
critical core mass, Menace ~ several M^, before the disk gas depletion. There are at least two factors 
which limit the gas accretion rate. Planetary atmosphere is no longer in hydrostatic equilibrium 
and gas accretion starts when Mc exceeds Mc^hydro — lO(Mc/lO~®M0yr~-^)''-^^M0 (Ikoma et al. 
2000). Pollack et al. (1996) suggested that M^. hydro > lOiW® becausc Mc may be maintained at a 
modest value, for a few Myr. But, cores emerge with neighbors of comparable masses. Together 
they stir up the random motion of the residual planetesimals, widen their feeding zone, and reduce 
their M^. In addition, the combined effects of gravitational scattering by the cores and gas drag 
actually lead to a clearing of residual planetesimals in the feeding zones of the cores (Tanaka & Ida 
1997; Rafikov 2003) and a reduction in 

hydro- Thus, the ousct of efficient gas accretion would 

occur earlier than previously estimated. 

With a negligible Mc, the gas accretion rate is determined by the planet's total (core plus 
envelope) mass (Mp) such that dMp/dt = Mp/rxH where We adopt 



with b = 10 and c = 3, which is a fitting formula for Pollack et al. (1996) 's results with negligible 
planetesimal accretion, tkh for Mp I-IOM0 can be longer if modest planetesimal accretion 
induced by gas accretion is assumed [e.g., "phase 2" by Pollack et al. (1996)]. But, the planetesimal 
accretion may be significantly reduced as discussed above. In this case, phase 2 cannot be sustained. 
On the other hand, if the core's migration is very fast, the growth of a core does not decline until 
it acquires mass ~ Mc^noiso- In this case, phase 2 does not exist either (Alibert et al. 2004). So, in 
the nominal case, wc adopt (6, c) = (10,3), but (6, c) = (9,3) and (11,3.5) are also examined (the 
latter corresponds to the gas accretion including phase 2). 

The gas accretion rate increases rapidly with Mp. Such runaway gas accretion is terminated 
when the residual gas is depleted globally on timescale rjep or it is severely depleted locally in the 
vicinity of their orbits due to gap opening, which is assumed to occur if their Hill radius exceeds 1.5 
times the disk scale height (Paper I). When Mp becomes a significant fraction of Mj (Jupiter mass), 
planets' tidal torque on the disk becomes larger than the viscous torque of the disk gas, resulting 
in a partial gap opening, e.g., (Lubow et al. 1999). It reduces gas accretion rate. Bondi accretion 
limit also becomes important in the reduction. However, the elongated accretion timescale may be 
still shorter than tkh for Mp ~ several M^ , so that the limited accretion at high Mp may not affect 




(5) 




(6) 
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the mass distribution of gas giants significantly. The critical mass -Menace for actual formation of gas 

giants corresponds to the value of Mp for which tkh ~ '^dcp- Equation (6) shows Menace ~ several 
M0. From eqs. (4) and (3), we find that, for sufficiently large /d, it is possible for Mc(Tdep) > -Menace 
and Mc^iso > Afc,acc so that gas giant can form readily. 

After the gap opening, the planets undergo "type-II" migration, in contrast to the cores' "type- 
I" migration without a gap (Ward 1997). In order to compute the protoplanets' migration, we adopt 
the a-prescription for effective viscosity (Shakura & Sunyaev 1973) with a uniform a = 10~^ for all 
the disks such that their viscous diffusion time scale (r,^) at a ~ lOAU is comparable to r^ep = 1-10 
Myrs. The planets migrate (Lin &; Papaloizou 1985, 1993) with a viscously evolving disk on time 
scales (Paper I) 




The planets' migration is terminated either when the disk is severely depleted or when they reach 
the vicinity of their host stars. We artificially halt the planetary migration at 0.04AU. 

Using these models, we carry out Monte Carlo simulations. For simplicity, we generate a set 
of initial a's of the protoplanets, T^ep a^nd with uniform distributions in log scale in the ranges 
of 0.1-lOOAU, 10^-10^ yrs, and O.7-1.4M0, respectively. The assumed a distribution corresponds 

to orbital separations Aa that is proportional to a, which may be the simplest choice. More refined 
choice of Aa that takes into account interactions between planets and detailed formation process 
will be explored in a future paper. However, as long as Aa increases with a with a moderate power 
index, calculated distributions of gas giants may not change significantly. 

We discuss the distribution of /d and fg^ in details, because it is a key assumption of this 
paper. A large range of dust mass in a disk is inferred. Observationally, the mass transfer rate from 
disks onto classical T Tauri stars is ~ 10~^-10''^Mq yv"^ (Calvet et al. 2000) with inferred disk 
masses equivalent to 0.1 < /g,o ^ 10. Independent mm observations (Beckwith &; Sargent 1996) 
show a similar range of fd- Here we assume that log^Q fgfl has a Gaussian distribution centered at 
with the dispersion 1.0 and upper cut-off at fg^ = 30 (Paper I). 

Since in the MMN, /d = /g,o = 1) the independently determined mean values (/d) ^ (/g,o) ^ 1 
imply that the average value of metallicity in the disk (Zdisk) is ^ Zq (the Solar composition), 
where [.^disk = logio(/d//g,o)]- In Paper I, /d = /g,o was assumed. Here we extend our calculation 
to a range of .^disk in order to account for the most pronounced correlation between the percentage 
of stars with planets detectable with the current Doppler survey and metallicity of their host stars 
(Z^) (Fischer Sz Valenti 2003; Santos et al. 2004), assuming ~ .^disk- From this inference, we 
explore the possibility that (/d) depends on as (/d) = lO^disk--^© ^ iqZ,-Zq^ while (/g,o) is ~ 1, 
independent of Z*, since gas disk formation/evolution processes may not be significantly affected 
by .^disk (~ Z^). We assume /d = /g,o x 10^* ~^®. This estimate can be verified with additional 
direct measurement of both gas and dust masses in protostellar disks in the near future. All stellar 
material is accreted through protostellar disks where differential retention of gas and dust would 
affect the stellar metallicity. These key assumptions are supported by the remarkable chemical 
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homogeneity among young stars in the Pleiades cluster (Wilden et al. 2002) and mature binary 
systems (Desidera et al. 2004). 

3. Numerical Results 

The simulated Mp-a distributions with AZ^ = Z^, — Zq = are shown in Figure la. This 
is the nominal case where m = 10^^ g, Mc,iso and tkr = W^^ {Mp / M(q)~'^ yrs are adopted. For 
comparison, the observed distributions are also plotted in Figure lb. The observational limits in 
the current Doppler survey, Mp > 100(a/lAU)^/2jv^^ a,nd a < 3AU (around F, G, K dwarfs within 
~ 50 pc) are marked by dotted lines. A deficit of planets with an asymptotic mass Mp = lO-lOOM® 
and a < 3AU (within a "planet desert") is apparent in Fig. la. This sparsely populated region 
is primarily caused by the runaway gas accretion onto cores with Mp > several M® (Paper I). 
This region divides the terrestrial (rocky) planet, gas and ice giant domains in the planets' mass- 
period distribution. Planetary migration generally sharpens the boundaries of the domains (Paper 
I; Udry et al. (2003)). Since migration is slower for larger Mp (Eq. [7]), relatively massive planets 
are less likely to migrate to inner regions, resulting in a deficit of planets with Mp > 1O^M0 at 
a < 0.2AU (Fig. la; also see discussion by Udry et al. (2003)). If larger a is adopted, becomes 
short compared with r^cp, so that even relatively massive planets can migrate significantly, resulting 
in less clear deficit. We found that the two planct-dcplctcd regions exist for other AZ* as well, 
because the "planet desert" and the other are determined by the core accretion model for gas giants 
and Tdep/Tj/, respectively. The observed distributions may exhibit the above two deficits and are 
consistent with our theoretical models with r^ep ~ Tj/ at ~ 10 AU. Adopting the same condition, 
the simulated a-distribution of the gas giant planets can match with the observation [Figure Ic; 
also see Armitage et al. (2002) and Trilling et al. (2002)]. 

The 0.04 AU cut-off in the planets' a distribution in our results may also be less abrupt and 
more consistent with the observed distributions if we adopt more relaxed inner boundary conditions. 
Many of giant planets may have migrated toward the vicinity of their host stars (Lin et al. 1996). 
The theoretically determined ratio of the giant planets with a < 0.06AU (hot Jupiters) to those 
between 0.2-2AU is ~ 1-10 for stars with Zq, which is order of magnitude larger than its 

observed value (~ 0.2). This discrepancy suggests that more than 90% of the hot Jupiters which 
migrated to the stellar vicinity are either consumed (Sandquist ct al. 1998) or tidally disrupted (Gu 
et al. 2003) by their host stars. Provided these disruptive events occur before the host stars become 
main scqTicncc stars with relatively shallow convection zones, their apparent chemical homogeneity 
may be preserved to match the observation of the Pleiades cluster (Wilden et al. 2002). As the 
migrating protoplanets sweep planetesimals along their paths, this process may self-regulate the 
amount of residual heavy elements in the disk. 

Figures Ic and d show that the normalized slopes of the a and Mp distributions are independent 
of AZ^, because they are determined by Tdep/Ti, and the core accretion model for gas giants, 
respectively. In these figures, we do not include planets stalled at 0.04AU, since a large fraction of 
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them may plunge into their host stars. 

Although the normalized distributions are independent of AZ*, we found that formation and 
retention rates of gas giant planets rapidly increase with AZ*. In order to account for the IS.Z^- 
dependence in the Mp and a distribution, we first describe how planet formation depends on [for 
details, see Paper I and Kokubo & Ida (2002)]. In a disk similar to the MMN (/d = /g^o = 1)) for 
a < Cice, the core growth is limited to Mc^iso < where is the Earth mass. Although growth 
beyond isolation (up to ~ M^) is possible for terrestrial planets, it is attainable only after severe 
gas depletion (Kominami &; Ida 2002). At these regions, both Mc^iso and Mc^noiso are < Menace such 
that the cores cannot acquire a massive gaseous envelope. At the outer regions where a !» Oice, 
although Mc,iso > M:,acc; the timescale required for Mc to reach Mc,acc is T^cp. The cores evolve 
into ice giants, similar to Uranus and Neptune. In the disk similar to th MMN, gas giants can 
form only slightly outside the ice boundary, because the cores may acquire -M^acc there within a 
few Myr and enter a rapid gas accretion phase prior to the global disk depletion. Gas giants can 
form inside Oice in some massive disks (with /d > 5). Rocky cores are formed with Mc^iso ^ in 
the inner regions and ice giants are formed with Mc(rdep) ^ in the outer regions of low-mass 
disks with /d 0.5, without formation of gas giants. 

Around metal-deficient stars, formation efficiency of gas giants is low because disks with suf- 
ficiently large /d are only at the tail end of the distribution. Around the metal-rich stars, the 
efficiency is high because both dM^/dt and Mc^iso increase rapidly with /d (eqs. 3 and 4) and 
-^c,iso ^ -^c,acc even interior to aicc- These account for the AZ*-dependence. In order to directly 
compare with the observed frequency of extrasolar planets, we plot in Figures 2 the theoretical de- 
termination for fraction (r/j) of stars which bear giant planets detectable with the current Doppler 
survey, as a function of their metallicity. Here, the distribution of /d is the same as that in Figures 
1. For each /d, a is selected as log(aj+i/aj) = 0.2 (j = 1,2, ...). We do not include planets stalled 
at 0.04AU in the evaluation of ryj. 

These results indicate ?7j increases linearly with AZ*, which is consistent with observational 
results of Fischer & Valenti (2003) (open circles in Figures 2). As shown in Figure 2a, ?7j does 
not depend on the choice of tkh [(ft, c) = (9, 3), (10, 3) or (11, 3.5); trh = 10^(Mp/Me)"''] nor that 
of asymptotic core mass [min(Mc,iso, Mc,noiso) or Mc^noiso]- The reasons are as follows. As long as 
Tkh has strong dependence on Mp, even significant change in h results in relatively small change 
of Menace, the core mass that enables rapid gas accretion prior to disk gas depletion. Since in our 
simulations, most gas giant planets are formed with such parameter ranges where /d and a are 
relatively large and Mc^iso ^ Mc^noiso, the distributions hardly change even if Mc^noiso is adopted in 
entire /d and a (although it affects the distributions of terrestrial planets [Paper I]). 

Although the absolute values of 77j are higher /lower for faster /slower core accretion, the depen- 
dence on AZ* is similar (Figure 2b). If disk depletion is caused by viscous diffusion, the assumption 
that Tdep ~ at ~ 10 AU is reasonable. However, the observation of dense clusters suggest that 
gas depletion in protostellar disks is primarily driven by photo-evaporation (Johnstone et al. 1998) 
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in these clusters, where a few massive stars contribute to nearly the entire UV ionizing flux. In 
order to consider the possibility of rapid gas depletion in a dense cluster, we carried out the addi- 
tional simulations with r^cp = 10^-10^ yrs. In this case, formation rates of giant planets are lower, 
resulting in smaller rjj (Figure 2b). But, the dependence on AZ* does not change also in this case. 
This general agreement is a generic feature of the core accretion scenario and the assumption for 
/d and /g,o- It may not arise naturally from the gravitational instability scenario (Boss 1997). 

4. Discussion 

Through Monte Carlo simulation using the deterministic planet formation model developed 
in Paper I, we find that the formation probability of gas giant planets increases rapidly with the 
metallicity of their host stars. Our model quantitatively accounts for observational correlation 
between stellar metallicity and the fraction of stars harboring giant planets that are detectable 
with the current Doppler survey. Our model also explains the two observationally suggested planet- 
depleted regions in Mp-a distributions of extrasolar planets: the deficit in intermediate mass Mp = 
10-lOOMe and a < 3AU, and that in large mass Mp > lO^M® and a < 0.2AU. 

Our model can be tested by more detailed observation. We show that the normalized slope 
of the a and Mp distributions of gas giants and the planet-depleted regions are independent of 

AZ*, because these are regulated by t^jcp/ti/ and the core accretion model. The lower limit mass 
lim of hot Jupiters depend on AZ^ as Mp Um oc 10 (3/2) AZ* (Paper I). The outer boundary 
radius of the distribution of gas giants Cout ^ ^q^"'^'^/^''^^^* (Paper I), although interactions among 
gas giants, which are not taken into account in the present paper, may make the boundary less 
clear (Marzari & Weidenschilling 2002). These predictions of AZ*-dependence can be tested by 
statistics of extrasolar planets with increased detection number in a few years. 

Our results provide supporting evidences for the core-accretion scenario and they suggest that 
the detection of gas giant planets may be used to infer the presence of rocky planets around the 
same host stars. "Habitable planets" which can retain liquid water on their surface (Kasting et al. 
1993) may have O.IM® < Mp < lOM® and a - 0.8-1.5 x (MjMo)^ AU. The results in Figure 1 
indicate that the frequency of " habitable planets" is comparable to or more than that of detectable 
giant planets (at present, more than 100 planets have been detected around F, G, K dwarfs within 
50 pc), although some fraction of the habitable planets may be destabilized by migrating giant 
planets. 
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FIGURE LEGENDS 
Fig.l 

Theoretically predicted mass (Mp) and semi-major axis (a) distributions of extrasolar planets. 
For parameters, see text. In panel (a), the computed result with AZ* = is shown. Totally 20000 
planets are calculated. For comparison, observational data are plotted in panel (b). (Since orbital 
inclinations are unknown, minimum masses are plotted). The observational limits in the current 
Doppler survey are marked by dotted lines. One dimensional normalized distributions of observable 
planets is plotted as a function of a (panel c) and Mp (panel d). Large open circles are observed 
results. Small filled circles, pentagons, squares, triangles, and crosses express the theoretical results 
with AZ^ = 0.5, 0.25, 0, —0.25 and —0.5, respectively. Here 1000 observable planets are used in each 
AZ^ result. 

Fig.2 

The fraction (r/j) of stars which bears giant planets currently detectable with Doppler survey, as 
a function of their metallicity (AZ*). Large open circles are observational results (Fischer &i Valenti 
2003). Other small symbols with lines represent theoretical predictions. Filled circles express the 
nominal case where core accretion is given by Eq. (3) with m = 10^^ g, the asymptotic masses 
are min(Mc^iso, Mc^noiso)) and tkh = ^O^iMp/M^)^^ yrs with (6, c) = (10,3). In panel (a), crosses, 
filled triangles and squares express the cases with Mc^noiso {b,c) = (9,3), and (6, c) = (11,3.5), 
respectively. In panel (b), crosses, filled triangles and squares express the cases with 3 times faster 
core accretion rate, 3 times slower one, and rjep = 10^-10^ yrs, respectively. 
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